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We study synchronization in delay-coupled neural networks of heterogeneous nodes. It is well 
known that heterogeneities in the nodes hinder synchronization when becoming too large. We 
show that an adaptive tuning of the overall coupling strength can be used to counteract the effect 
of the heterogeneity. Our adaptive controller is demonstrated on ring networks of FitzHugh- 
Nagumo systems which are paradigmatic for excitable dynamics but can also - depending on 
the system parameters - exhibit self-sustained periodic firing. We show that the adaptively 
tuned time-delayed coupling enables synchronization even if parameter heterogeneities are so 
large that excitable nodes coexist with oscillatory ones. 


1. Introduction 


The ability to control nonlinear dynamical systems has brought up a wide interdisciplinary area of research 
that has evolved rapidly in the past decades |Scholl h Schuster, 20081. Besides the control of isolated 
systems, control of dynamics in spatiotemporal systems and on networks has recently gained much interest 
|Sun fe Ding , 2013; Flunkert et al ., 2010; Omelchenko et al., 2011; Postlethwaite et al., 2013; Bachmair fc 
Scholl, 2014| . Adaptive control schemes have emerged as a new type of control methods that enable control 
in situations where parameters are unknown or drift in time. They allow for automatically tuning the 
parameters to appropriate values and are, therefore, of particular interest for experiments and technological 
applications. Previously, they have been successfully applied in the control of network dynamics |Selivanov 


et al ., 2012 Lehnert et al. , 2014 . Here, we show that adaptive control methods can be used to counteract 


the effect of heterogeneous nodes in the synchronization of delay-coupled networks. 

Synchronization in neural networks has gained a lot of attention lately [Pikovsky et al. 


is involved in processes as diverse as learning and visual perception on the one hand ies| 
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Guzenko et al. , 2013; Lehnert et al. , 2014| . Further results on the control of network dynamics were obtained 


in |Lu et al. , 2015; Li et al. , 2014] . However, in realistic networks the nodes will always be characterized 
by some diversity meaning that the parameters of the different nodes are not identical but drawn from a 
distribution. It is well known that such heterogeneities in the nodes can hinder or prevent synchronization 
and that the coupling strength is a crucial parameter in this context |Strogatz, 2000; Sun et al., 20091. 
Moreover, in most cases perfect synchronization - in the sense that the state of all nodes is identical at all 
times - is unfeasible in the presence of heterogeneities. Although quite a number of papers are devoted to 


synchronization of heterogeneous networks, see Refs. |DeLellis et al. , 2015; Isidori et al ., 2014; Chen, 2014 


Ricci et al. , 2012] and references therein, only a few of them address adaptive control of synchronization 


|Guo fc Li , 2013; Fradkov & Junussov, 2013 . Moreover, the adaptation of the coupling strength was not 
considered previously in the case of heterogeneous networks. Here, we develop an adaptive controller which 
allows for reaching a synchronization state where the node’s trajectories are close though not identical. 

Our method is based on the speed-gradient (SG) method, which was previously used in the control of 
delay-coupled networks |Selivanov et al ., 2012; Guzenko et al ., 2013; Lehnert et al ., 2014] , however, not in 
the presence of node heterogeneities. In order to apply the SG method, we suggest a goal function which 
characterizes the quality of synchronization. Based on this measure an adaptive controller is developed 
which ensures synchronization even if the parameter heterogeneities become such large that some nodes 
- if uncoupled - undergo a Hopf bifurcation and behave distinctly different from the other nodes in the 
network. We demonstrate our algorithm on the FitzHugh-Nagumo (FHN) system |FitzHugh, 1961; Nagumo 


et al., 1962 , a generic model for neural dynamics. The advantage of our approach is that the SG method 


allows for using a simple goal function which does not depend upon the parameters to be controlled but 
only on the system state. 

Note that in contrast to artificial neural networks Refs. |Li et a/." , 2014; Lu et al ., 2015; Shi et al. 


2015 which are designed to solve advanced computational tasks, we here assume that the network, i.e., the 


local dynamics and the network topology, are given and only design the algorithm adapting the coupling 
strength. Furthermore, usually all nodes in artificial neural networks are identical, while we try to mimic 
realistic networks and therefore use heterogeneous network nodes. 

The paper is organized as follows: Section[2]is a recapitulation of the SG method, while Sec. [^introduces 
the model. Section [4] discusses two delay-coupled FHN systems: The possible bifurcation scenarios are 
investigated and the adaptive control algorithm is developed. In Sec. [5| the method is generalized to larger 
ring networks. Finally, we conclude with Sec. [6j 


2. Speed-Gradient Method 


In this section, we briefly review 
nonlinear dynamical system 


the speed-gradient (SG) method [Fradkov, 2007 . Consider a general 


x — F(x, g, t) (1) 

with state vector x G C n , input (control) variables g G C m , and nonlinear function F. Define a control goal 

Q(x(i),i)< A, (2) 

for t > £*, where Q(x,t) > 0 is a smooth scalar goal function and A is the desired level of precision. For 
example, if we want to force the trajectory of system 0 to follow the desired trajectory x*(£), we can use 
a goal function in the form Q(x(t)) = (x(£) — x*(£)) 2 . 

In order to design a control algorithm, the scalar function Q = cj(x, g, t) is calculated, that is, the 
speed (rate) at which Q(x(t),t) is changing along the trajectories of Eq. Q: 

<^(x, g, t) = - + [V x Q(x, i)] T F(x, g, t). (3) 

Then the gradient of cj(x, g, £) with respect to the input variables is evaluated as 

V g cj(x,g ,t) = V g [V x Q(x, f)] T F(x, g, t). 


( 4 ) 
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Finally, we obtain the control function g from 

g(*) = g°-'i/’(x,g ,t), 


( 5 ) 


where the vector function ^(x, g,t) = 7 V g w(x,g,t) with some adaptation gain 7 > 0, and g° = const is 
an initial (reference) control value (often g° = 0 is assumed). The algorithm © is called speed-gradient 
(SG) algorithm in finite form since it suggests to change g proportionally to the gradient of the speed of 
changing Q. For the speed-gradient (SG) in its differential form see Ref. |Fradkov, 2007 . 

Several analytic conditions exist guaranteeing that the control goal (|2j) can be achieved in system 0 
and ©. The main condition is the existence of a constant value of the parameter g*, ensuring attainability 
of the goal in the system dx./dt — F(x, g*, t). Details can be found in the control-related literature |Fradkov[ 
1979| |Shiriaev fe Fradkov] 2000 | . 


The idea of this algorithm is the following: The term — V g cj(x, g, t) points to the direction in which 
the value of Q decreases with the highest speed. Therefore, if one forces the control signal to ’’follow” this 
direction, the value of Q will decrease and finally be negative. When Q < 0, then Q will decrease and, 
eventually, will tend to zero. 


3. Model equation 

The local dynamics of each node in the network is modeled by the FitzHugh-Nagumo (FHN) differential 
equations |FitzHugh , 1961; Nagumo et al , 1962| . The FHN model is paradigmatic for excitable dynamics 
close to a Hopf bifurcation jLindner et al , 2004| , which is not only characteristic for neurons but also occurs 
in the context of other systems ranging from electronic circuits | Heinrich et al , 2010 ] to cardiovascular 
tissues and the climate system |Murray , 1993; Izhikevich, 2000| . Each node of the network is described as 
follows: 


N 


em = m - 


Vi — Ui T cq 


- Vi 


+ C Yi Gi i K - r ) - U i (*)] > 

3 = 1 

i = 1,..., TV, 


( 6 ) 


where U{ and V{ denote the activator and inhibitor variable of the nodes i — 1,..., TV, respectively, r is 
the delay, i.e., the time the signal needs to propagate between node i and j (here we will use r = 1.5). e 
is a time-scale parameter and typically small (here we will use e = 0.1), i.e., U{ is a fast variable, while 
Vi changes slowly. The coupling matrix G = {Gij} defines which nodes are connected to each other. We 
construct the matrix G by setting the entry Gij equal to 1 (or 0) if the jth node couples (or does not 
couple) into the ith node. After repeating this procedure for all entries Gij, we normalize each row to unity. 
The overall coupling strength is given by C. 

In the uncoupled system ((7 = 0), ai is a threshold parameter: For a i > 1 the zth node of the system 
is excitable, while for cq < 1 it exhibits self-sustained periodic firing. This is due to a supercritical Hopf 
bifurcation at cq = 1 with a locally stable equilibrium point for a i > 1 and a stable limit cycle for cq < 1. In 
previous publications, networks of homogeneous FHN systems were considered, i.e., a\ — — ... = ajsi = a 


Brandstetter et al 2010 

; Scholl et al , 2009 

Hovel et al 

, 2010; Lehnert et al/. 2011 

Panchuk et al 

2013 

Cakan et al 

, 2014 . In particular, it was shown that for excitable systems, i.e., a > 1 and coupling matrices 


with positive entries zero-lag synchronization is always a stable solution independently of the coupling 
strength and delay time (as long as both are large enough to induce any spiking at all). 


Here, we investigate the case of heterogeneous nodes. In this case, perfect synchronization, i.e., 
(^ 1 ,^ 1 ) = ... = (un,vn) = (u s ,v s ), is no longer a solution of Eq. 0 which can easily be seen by 
plugging (ui,vi) = ... = (u N ,v N ) = (u s ,v a ) into Eq. 0. The node dynamics is then described by 

£U S = u s - y - t) s + C2][ii s (i - r) - u s (t)}, ^ 

v s = u s + a,i , i — 1,..., N, 


which is obviously not independent of i. This means that a perfectly synchronous solution does not exist 
in system 0 because the prerequisite for the existence of such a solution is that each node receives the 
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same input if all nodes are in synchrony. However, solutions close to the synchronous solution might exist 
where the nodes spike at the same (or almost the same) time but with slightly different amplitudes. As we 
show, these solutions can be reached and stabilized by an adaptive tuning of the coupling strength. 

4. Two delay-coupled FitzHugh-Nagumo systems 

This Section studies the most basic network motif consisting of two coupled systems without self-feedback. 
Before deriving the adaptive controller, we perform a linear stability analysis of the equilibrium point to 
get insight in the possible bifurcations. 


4.1. Linear stability of the equilibrium point 

The linear stability analysis follows the approach suggested in Ref. |Dahlem et al, 2009; Scholl et aZ., 2009 . 
Consider two coupled FHN-systems with heterogeneous threshold parameters and bidirectional coupling 


£Ui = Ui - Y - V 1 + C[u 2 (t - t) - Ui(i)], 


v\ = u\ + ai, 
u 3 


( 8 ) 


£U 2 = U 2 
V2 = U2 P CL2- 


Y~ V2 + C i U ^( t ~ T )~ «2 (t)], 



/a -1 0 0 \ 


/o oco\ 

1 

£ 0 0 0 

r 1 

0 0 0 0 

£ 

0 0 ( 2-1 

5x(t) + - 

£ 

CO 0 0 


^0 0 e 0 ) 


^0 0 0 0/ 


^At- 


where ^ = 1 — a? — C. The ansatz 

5x(t) = e A ”q, 

where q is a time-independent vector, leads to the characteristic equation for the eigenvalues A, 

(1 — £iA + sA^)(l — £2 A + s\ 2 ) — (XCe Ar )^ = 0 . 


The unique equilibrium point of the system ([8]) is given by x* = (u\, v*, > where u\ — — ai, 

U 2 — —CL 2 , v* = — a\ + a\/ 3 + C{a\ — 0 , 2 ) and ^ — ~ a 2 + ^ 2/3 + C(d 2 — a\). Linearizing Eq. Q around 
the equilibrium point x* by setting x(£) = [iq(£), 'iq(t), U 2 (t), V 2 (t)] T = x* + 5x(t), we obtain 

5x(t-r), (9) 

( 10 ) 
( 11 ) 

A necessary condition for the occurrence of a Hopf bifurcation is that the eigenvalue A is imaginary. 
We, therefore, substitute the ansatz A = icj, uo G R, into Eq. © to find the threshold parameters a\ and 
d 2 for which a Hopf bifurcation can take place. Then, separating Eq. © into real and imaginary parts 
yields 

(1 - £UJ 2 ) 2 - = -lo 2 C 2 cos(2wt), 

w(l — ew 2 )(^i + £ 2 ) = — tz 2 C 2 sin(2o;r). 

Since £<1, terms of the order of e can neglected, i.e., we consider e — 0 in the following. Squaring 
and adding the subequations of Eq. ( [I2| ) we obtain 

<^ 2 (Ci + C 2) 2 + (1 ~ Ci^2<^ 2 ) 2 = c/C 4 , (13) 

which can be rearranged to 

(^ 2 2 - ^ V + (Ci 2 + &W + 1 = 0. (14) 

Equation (14) is biquadratic, i.e. quadratic in 2 ; = uj 2 . Therefore, we can use Vieta’s formulas to analyze 


whether it has non-negative roots. According to Vieta the following holds 

tf + tl 1 


Zi + *2 = 




ziz 2 = 


ep 2 -c^ 


(15) 
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where z\ and z 2 are roots of the quadratic equation (14). For > C 4 , z\ < 0 and Z 2 < 0 follows meaning 
that Eq. (14) has no real-valued solution cj, and, thus, no Hopf bifurcation will take place. Taking the 
square root of this inequality and resubstituting ^ = 1 — a? — C yields 


\(l-C-a 2 1 )(l-C-a 2 2 )\>C 2 . 


(16) 



(a) (b) 

Fig. 1. Hopf bifurcation of the equilibrium point of two delay-coupled FHN systems (Eq. §) for (a) small coupling strength 
C, i.e., C — 0.3, and (b) large coupling strength C — 5. (Red) shading marks parameter values for which inequality (16) is 
fulfilled, i.e., a Hopf bifurcation is impossible and the equilibrium point is stable. 


This inequality defines the values of a\ and <22 where a Hopf bifurcation is impossible, i.e., the equilib¬ 
rium point is stable. Figure [l] shows this area in (red) shading for (a) small coupling strength, i.e., C < 1, 
and (b) large coupling strength, i.e., C > 1. As expected no Hopf bifurcation will occur for a\ > 1 and 
<22 > 1 since then both systems are in the excitable regime. Note that, however, oscillations can coexist 
with the stable equilibrium point |Lehnert et al. , 2011 . 


4.2. Adaptive control of two coupled FHN-systems 

We now want to apply the SG method to system ®> with the goal to synchronize the two heterogeneous 
nodes. As discussed above perfect synchronization in the form (ui,vi) = (^ 2 ,^ 2 ) is not attainable in this 
case but the two systems will follow slightly different trajectories in the synchronized case. We, therefore, 
use as a goal function 


1 




(17) 


u 2 (t) + 01 - a 2 ) 2 . 

The choice ensures that the system follows trajectories for which 

ui(t) - u 2 (t) « —<21 + a 2 , (18a) 

v\ (t) - V 2 (t) « c (18b) 

holds for t > £*, where c is a constant. Approximati on (|18a ) directly follows from the chosen goal function 


( [T7| ); approximation (18b) is obtained by plugging ( 18a| ) into Eq. ([§]). Thus, the goal function ( pT| yields 
synchronization with a shift in the values of the inhibitors and activators of the two nodes. 

From Eq. © with g = C, system Q, goal function <HZt and ^(x, (7, t) = 7 Vc^(x, (7, t) an adaptive 
law is straightforwardly derived: 


C(t) = Co+ 1 [«1 (t) 
£ 


U 2 (t) + a\ - a 2 \ [ui(t) - U 2 (t) +ui(t — r) — u 2 {t - r)] 


(19) 


where 7 > 0 is the gain and Co is the initial value of the control parameter. The appropriate value of 7 
has to be determined by numerical simulations. Note that a similar approach has been used to tune the 
coupling strength in a network of Rossler systems in Ref. (Guzenko et al. , 2013]. 

Close to the control goal, u\(t) — U2(t) ~ u\{t — r) — U2(t — r) ^ a\ — <22 holds. We, therefore, can 
simplify the adaptation law by substituting the delayed variables by their non-delayed versions and obtain 

C(t) = Co + — [u\(t) - u 2 (t ) + ai - a 2 ] [ui(t) - u 2 (t)\. 


( 20 ) 
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Fig. 2. Dynamics of two coupled FitzHugh-Nagumo systems according to Eq. j8| with constant coupling strength C. Gray 
(green) line marks node one, black (red) line marks node two. (a) and (b): time series of the activator and the inhibitor, 
respectively; (c): difference u\ — U2 between the activator values, and (d): phase space. Parameters: N = 2, e — 0.1, r — 1.5, 
ai = 1.1, a 2 — 0.7, C — 1. Initial conditions: iq(£) — Vi(t) = 0, i— 1,2, for t £ [— r, 0]. 


For constant coupling strength, i.e., 7 = 0, the two coupled FHN systems do not synchronize in-phase, 
but approach an anti-phase synchronized state: Figure [2] shows in panels (a) and (b) the time series of the 
activators and the inhibitors, respectively, and in panel (d) the phase portrait. Though the first node is in 
the excitable regime (ai = 1.1 > 1) both nodes oscillate due to the nonzero coupling strength C. However, 
they do not synchronize as can clearly be seen in panel (c) which depicts the difference u\ — U 2 between 
the activator values. Instead, they phase lock with a phase shift of approximately n which corresponds to 
an anti-phase synchronized state. 


We now adapt the coupling strength according to Eq. (20) in order to synchronize the two systems 
where the result is shown in Fig. [3| After a transient time of approximately 15 units of time the two systems 
reach the desired synchronized state (see the time series of the activator in Fig. ^a) and the difference 
between its values Fig. [3](c)). Thus, the control is successful. 

If the gain 7 is chosen too low the control fails: Figure [4] depicts the results of the adaptive control 
according to Eq. ( |20| ) for 7 = 0.05. Clearly, the control does not succeed in synchronizing the two systems 
(see the time series of the activators in Fig. Qa) and the difference between their values in (c)). 

The adaptive controller (19) ensures synchronization of the activators with a shift given by <22 — a\ 
(see Eq. (18a)). Furthermore, there is a finite, constant shift in the inhibitor values (see Eq. ( |18b| )). The 
shift in the inhibitor values can be reduced to a value close to zero if we control the coupling strength of 
each node separately. The two FHN systems are then described by 




EUi = Ui - y - Vi + Ci(i)[U( i+1 ) mod 2(^ “ t) ~ Ui(t )], 


( 21 ) 


Vi — Ui -\- cii 


i = 1,2 


where Ci(t ) describes the strength of the coupling to node i. From Eq. © with g = (Ci, C 2 ), system ( [ 21 ] ) , 
goal function < 0 * and ^(x, g ,£) = 7 V g cd(x, g,£) an adaptive law is straightforwardly derived: 

Ci(t ) C* + — [rq(t) 2 £p_|_ 1 ) mod 2 ^) T CLi ^(«+l) mod 2] [^i(^) ^(z+ 1 ) mod 2 ^ 5 ^ ^ ( 22 ) 

where 7 > 0 is the gain and C® is the initial value of control parameter. 

The results of the adaptation according to Eq. (22) are shown in Fig. [5j The two systems reach the 
desired synchronized state (see the time series of the activators in Fig. [5](a) and the difference between 
their values in Fig. [5^c)). Moreover, the difference between inhibitor values is close to zero (see the time 
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(e) 


(f) 


Fig. 3. Adaptive control of two coupled FitzHugh-Nagumo systems (Eq. <©)• Gray (green) line marks node one, black (red) 
line marks node two. (a) and (b): time series of the activator and the inhibitor, respectively; (c) and (d): differences u\ — U 2 
and v\ — V 2 between the activator and the inhibitor values, respectively; (e): phase space, and (f): time series of the coupling 
strength adapted according to Eq. (20). Parameters: 7 = 3, Cq — 0. Other parameters and initial conditions as in Fig [ 2 ] 




(a) (b) 




(c) (d) 

Fig. 4. As in Fig. [ 3 ] but with smaller adaptation gain 7 = 0.025. Gray (green) line marks node one, black (red) line marks 
node two. (a) and (b): time series of the activator and the inhibitor, respectively; (c) and (d): differences u\ — U 2 and v\ — V 2 
between the activator and the inhibitor values, respectively. 
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t 



(a) 


(b) 



Fig. 5. Adaptive control of two coupled FitzHugh-Nagumo systems (Eq. (21)) where the coupling strength Ci of each node 
is adapted separately according to Eq. (22). Gray (green) line marks node one, black (red) line marks node two. (a) and (b): 
time series of the activator and the inhibitor, respectively; (c) and (d): differences u\ — U 2 and v\ — V 2 between the activator 
and the inhibitor values, respectively. Parameters: 7 = 2, C® — 0, i = 1, 2. Other parameters and initial conditions as in Fig [2] 


series of the inhibitors in Fig. [5](b) and the difference between their values in Fig. [ 5 pl)). Thus, the control 
is successful. 


5. Adaptive Synchronization in ring networks 

We now want to apply our method to larger networks. To this end, we consider a ring network of N nodes 
where the coupling matrix G has the following form 


fOl 0 ■ 
00 1 • 


G = 


00 

Vio 


0 


°\ 

0 


. 1 
•0 / 


(23) 


Our approach aims to synchronize two nodes of the ring; in this case the other nodes will follow these 
two nodes and synchronize as well. The same idea has been used in the control of wave motion in a chain 
of pendulums | Fradkov & Andrievsky] [2008 . However, there is a limitation to this approach: The nodes 
with the highest threshold and the lowest one have to be neighbors. 

Let us assume that k is the node with the highest threshold and l is the one with the lowest threshold, 


i.e, 


dk = max 07 , 


ai = mm 07 , 
z=l,...,iV 


and that a & and a/ are neighbors, i.e., k = (1 + 1) mod N or k = (l — 1) mod N. 

We now use the adaptation law (20) to synchronize these two nodes. With nodes k and l instead of 


nodes 1 and 2, Eq. (20) reads 


7 


C{t) = Co + - [' U k (t ) - Ui(t ) +a k - at} [ u k (t ) - u t (t )] , 
e 


(24) 


where 7 is the gain and Co is an initial value of the coupling strength. As before, we aim to achieve the 
control goal <§• Similar to the case of two coupled FHN systems, the control ( |24| ) ensures synchronization 
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(a) (b) 



(e) (f) 


U 



t 


(g) 



t 


(h) 

Fig. 6. Adaptive control of synchronization of a ring of ten FitzHugh-Nagumo systems (according to Eq. © with coupling 
matrix (23)). Gray (green) line marks node one, black (red) line marks node two. (a) and (b): time series of the activator 
and the inhibitor, respectively; (c) and (d): differences between the activator and the inhibitor values, respectively; (e): phase 
space, (f): time series of the coupling strength adapted according to Eq. (24); (g) and (h): time series of the activators and 
the inhibitors of all nodes, respectively. Parameters: N = 10, e — 0.1, r = 1.5, 7 = 100, Co = 0. Threshold parameters ai are 
chosen randomly from [0.8,1.1], while the highest threshold, here 01, equals 1.1 and the lowest threshold, here <22, equals 0.8. 
Initial conditions: Ui(t) = 0, Vi(t) = 0, i — 1, ..., N, for t G [—r, 0]. 


in a nearly synchronized state characterized by 


Ui(t) — Uj(t ) ~ — di + aj, 
vi(t) - vj(t) « a, 


(25a) 

(25b) 
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(a) 
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(b) 

Fig. 7. As in Fig.[6]but with a bounded coupling strength \C\ ^ 5. (a) and (b): time series of the activators and the inhibitors 
of all nodes, respectively. 




(b) 

Fig. 8. Adaptive control of synchronization of ten FitzHugh-Nagumo systems where the coupling strength of each node is 
adapted separately according to Eq. (29). (a) and (b): time series of the activator and the inhibitor of all nodes, respectively. 
Parameters: Cq = 0, i — 1,..., IV, 7 = 1. Other parameters and initial conditions as in Fig. ^ 


for t ^ £*, where <7 is constant and i, j = 1 ,..., N. 

Figure [ 6 ] presents the results of a simulation of ten FitzHugh-Nagumo systems coupled in a ring where 
the adaption law ( |24| ) is applied with a\ — a & being the node with the highest threshold, and a 2 = a\ being 
the node with the lowest threshold. In Fig. [ 6 ](a) and[ 6 ](c), it is shown that activators u\ and U 2 synchronize. 
Figure [ 6 ^e) shows the phase space of the two nodes. In Fig. ^g) it can be seen that not only node one 
and two are synchronized but that all nodes follow these two nodes and synchronize after a transient time. 
Thus, the control goal is achieved. 

Figure [ 6 ^f) shows that there are several jumps in the control C during the transient time and that 
C temporarily reaches quite high values. To avoid this, we can limit the coupling strength C. Figure [ 7 ] 
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presents the results for a bounded coupling strength, i.e., \C\ ^ 5. The control is successful; however, the 
transient time is prolonged compared to the unbounded control shown in Fig. [ 6 ^g). 

So far we considered the case that the nodes with the highest and lowest coupling strength are neigh¬ 
bors. If this is not fulfilled, the adaptation of the overall coupling strength does not yield synchronization. 
However, if we control the coupling strength of each node separately the control goal can be reached. The 
ring network is then described by 

,3 


u j 


EUi = Ui - y - Vi + Ci(t)[u( i+1 ) mod N (t-r) - Ui(t)}, 

Vi = Ui + ai, i = l,...,N, 


(26) 


where Ci(t) describes the strength of the coupling to node i. As in the case of two nodes, from Eq. ©> with 
g = (Ci,..., Cw), goal function 


1 N 2 
Q(x(t), t) — 'y ^ mod n(P) T &i &(%- (-1) mod N~\ i 


2=1 


and 


(27) 


(28) 


^(x,g ,t) = qV g w(x, g, t), 

we derive the following adaption law 

Ciit) =C° + J [uiit) - «(i +1) mod N (t - r)] (29) 

X ^ 2 rq(t) ^(2— 1 ) mod n(P) ^(2+1) mod n(P) T ^(i— 1 ) mod n(P) ^(2+1) mod iv(^)] 5 ^ • • • ? 

where 7 is the gain and C? is the initial value of C\. 

Figure [ 8 ] presents the results of a simulation of Eq. (26) with adaption law ( [29] ) . The control goal is 


achieved and synchronization takes place. This method requires a lower gain than algorithm (24), and, 
moreover, it ensures synchronization of the inhibitors with a negligible shift (see the time series of the 
inhibitors in Fig. | 8 ^b)). 


6. Conclusion 

We have proposed a novel adaptive method for controlling synchrony in heterogeneous networks. It is well 
known that networks with heterogeneous nodes are much less likely to synchronize than networks of iden¬ 
tical nodes. Furthermore, synchrony will take place in a state where the trajectories of the different nodes 
are not identical but small deviations can be observed. We have suggested a goal function to characterize 
this type of synchrony. Based on this goal function and the speed-gradient (SG) method, we have derived 
an adaptive controller which tunes the overall coupling strength such that synchrony is stable despite the 
node heterogeneities. 

We have demonstrated our method on networks of FitzHugh-Nagumo systems, a neural model which 
is considered to be generic for excitable systems close to a Hopf bifurcation. Before applying the adaptive 
control, we have studied the simple motif of two delay-coupled, heterogeneous nodes and have given analytic 
conditions for the occurrence of the Hopf bifurcations. We have then applied our adaptive method and 
discussed its dependencies on the node and control parameters. It has been shown that our method enables 
synchronization even if the node parameters are chosen such diverse that one of the systems would exhibit 
self-sustained oscillations without coupling, while the other one would remain in a stable equilibrium point, 
i.e., one of the uncoupled systems is above, and the other is below the Hopf bifurcation. Furthermore, we 
have generalized our method to larger networks and applied it to ring networks. As an alternative and 
complement to adapting the overall coupling strength, we have suggested adapting the coupling strength 
of each node separately. This allows for adaptively controlling the network in situations where the node 
with the highest threshold is not a direct neighbor of the node with the lowest threshold, in contrast to 
the restriction imposed by tuning only the overall coupling strength. 
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Here, we have considered static delays. In principle, it would be possible to extend the method to time 


Tang 

& 

-to 

QJ 

et al. 

2008; ( 


2008 


2008; Li et al. , 2014|. In general networks time varying delay has been studied in Refs. |Yu| 


Gjurchinovski et a/., 2014]. However, to our knowledge no adaptive control methods have been 


applied in networks with time varying delay and heterogeneous nodes. 

Given the paradigmatic nature of the FitzHugh-Nagumo system, we expect our method to be applicable 
in a wide range of excitable systems. Furthermore, the application of the SG method to the control of 
networks with heterogeneous nodes suggests that other adaptive controllers that are based on the SG 
method (see, for example, Refs. [Selivanov et a/. , 2012; Scholl et a/., 2012; Guzenko et a/., 2013; Lehnert 
et al. , |2014| ) are also robust towards heterogeneities. 
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